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Abstract 

§ Self-assembly is the mechanism that controls the formation of well defined structures from dis- 

ordered pre-existing parts. Despite the importance of self-assembly as a manufacturing method 

g 

I and the increasingly large number of experimental realizations of complex self-assembled nano 

C aggregates, theoretical predictions are lagging behind. Here we show that for a non-trivial self- 

I assembly phenomenon, originating branched loop-less clusters, it is possible to derive a fully pre- 

y—i dictive parameter-free theory of equilibrium self-assembly by combining the Wertheim theory for 

> 

associating liquids with the Flory-Stockmayer approach for chemical gelation. 
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Intermolecular self-assembly is the ability of molecules to form supramolecular assem- 
blies PP as well as a manufacturing method used to construct aggregate at the nano or 
micro scale, by proper design of the constituent molecules. In the self-assembly bottom- 
up paradigm, the final (desired) structure is 'encoded' in the shape and properties of the 
designed building blocks. Realizations of complex self-assembled nano aggregates [21 [3l S] 
have been guided by intuition and sophisticated experimental techniques. A full compre- 
hension of the self-assembly process requires the ability to predict the structures (and their 
relative abundance) which will be observed in equilibrium, as a function of T and density 
p, starting from the knowledge of the inter-particle interaction potential. Such a request 
is very much akin to the one that has guided the development of the physics of liquids in 
the last decades. Differently from the simple liquid case, self-assembly is characterized by a 
very strong inter-particle attraction (significantly larger than the thermal energy ksT) and 
by the fact that the interaction geometry is far from being spherical. The leading "bond- 
ing" interaction may indeed be localized in a specific part of the particle surface (patchy 
interactions [5j), it may be active only in the presence of a specific complementary group 
(lock-and-key interactions, very often encountered in biological self-assembly P [71 El E]) 
or it may be strongly dependent on the particle orientation [10]. The presence of strong 
and patchy interactions poses significant challenges to a parameter-free description of self- 
assembly. Only equilibrium chain polymerization, the simplest self-assembly process which 
takes place when bi-functional particles self-assemble into chains of variable length, can be 
considered to be sufficiently estabhshed [IH 121 ESI El US US] ■ 

In this article we show that for systems with a small average functionality (but larger 
than two) it is possible to provide a parameter-free full description of the self-assembly 
process. We study theoretically and numerically one of the simplest, but not trivial self- 
assembly process, namely a binary mixture of particles with two and three attractive sites. 
The presence of three (or more) -functional particles — which act as branching points in 
the self-assembled clusters — introduces two important phenomena which are missing in 
equilibrium chain polymerization: a percolation transition, where a spanning cluster appears, 
and a region of thermodynamic instability, the analog of a gas-liquid phase separation. More 
explicitly, we investigate a binary mixture composed by N2 = 5670 bi-functional particles 
and = 330 three-functional ones. The resulting average number of sticky spots per 
particle, i.e. the average functionality, is (M) = Particles are modeled as hard- 



2 



spheres of diameter cr, whose surface is decorated by two (or three) bonding sites at fixed 
locations. Sites on different particles interact via a square-well potential Vsw of depth uq 
and attraction range S — O.llQcr. More precisely, the interaction potential 1/(1,2) between 
particles 1 and 2 is 



where Vhs is the hard-sphere potential and ri2 and r^^^ respectively the vectors joining 
the particle-particle centers and the site-site (on different particles) locations; indicates 
the number of sites of particle i. Sites are located on the poles in the case of bi-functional 
particles and equidistant on the equator in the case of three-functional particles. With 
this site geometry, the smallest possible bonded ring requires six three-functional particles, 
while a bonded ring of only bi-functional units requires ~ 50 particles. The well width 5 is 
chosen to ensure that each site, due to steric effect, is engaged at most in one interaction. 
Distances are measured in units of a. Temperature is measured in units of the potential 
depth (i.e. Boltzmann constant ks = I)- In the studied model, bonding is properly defined: 
two particles are bonded when their pair interaction energy is -Uq. This means that the 
potential energy of the system is proportional to the number of bonds. The lowest energy 
state of the system (the ground state energy) coincides with configurations in which all 
bonds are formed, i.e Egg = —uq '^^'^'^^^^ . As a result, the bond probability, pb, can be 
precisely calculated as the ratio of the potential energy E and Egg. Pairs of bonded particles 
are assumed to belong to the same cluster. 

We have performed standard Monte Carlo (MC) Metropolis simulations at several T 
and p (more than 70 different state points). An MC step is defined as an attempted move 
per particle. A move is defined as a displacement of a randomly selected particle in each 
direction of a random quantity distributed uniformly between ±0.05 a and a rotation around 
a random axis of a random angle uniformly distributed between ±0.1 radiant. Equilibration 
was monitored via the evolution of the potential energy (a direct measure of the number 
of bonds in the system). Equilibration at the lowest T required up to 10^ MC steps (6 
months of simulation time). Indeed, at kBT/uo — 0.05, a bond breaks in average every 5 10^ 
attempts. 

The extremely long Monte Carlo simulations provide a numerically "exact" description 
of the equilibrium self-assembly process in this model. On cooling particles aggregate in 
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open larger and larger clusters which then coalesce into a spanning structure, eventually 
incorporating all particles. A pictorial representation of the structure of the system on 
varying T (both below, at and above percolation) is shown in Fig. [l] Clusters are composed 
by rather stiff chains of bi-functional particles cross-linked by the three-functional ones. 

To develop a parameter-free theoretical description of the cluster formation we combine 
Wertheim (W) theory pTl [T8l [19] for estimating the T and p dependence of pb with the Flory- 
Stockmayer (FS) predictions for the cluster size distribution in chemical gelation, providing 
a detailed and predictive theory of reversible self-assembly. 

The Wertheim theory, developed back in 1980 for describing the free energy of molecules 
with fixed valence (associating liquids), can be transferred to particles with a small number of 
patchy interacting sites. The main assumption in the theory is that molecules (or particles) 
cluster in open structures without closed bond loops. Such a condition, as we will show in 
the following, can be realized with patchy particles when the average functionality is small 
and the chains of bi-functional particles are not significantly flexible. The Wertheim theory 
predicts that can be calculated from the chemical equilibrium between two non-bonded 
sites forming a bonded pair. For the present model, such relation reads: 

(M)pA (2) 



where A ^ A. j ,Mm))...A.dn.. Here ...(12) .s the reference hard-sphere ra- 
dial distribution function, /(12) = e^^sw-{r/2)/'=s7' _ ^ jg ^j^g Mayer /-function between two 
arbitrary sites i and j, and (/(12))^^^^2 represents an angle average over all orientations of 
the two particles at fixed relative distance ri2 [20] . The comparison between the theoretical 
predictions and the "exact" numerical data for the T and p dependence of ph are shown in 
Fig. |2} Data show clearly that the theory is able to predict precisely pi, (or equivalently 
the system potential energy) in a wide T and p range. At low T, p;, — > 1, and the system 
approaches a fully bonded disordered (ground state) configuration. 

To derive information on the structure of the system and the connectivity of the aggre- 
gates we connect the W and the FS theories. Indeed, the hypothesis of absence of closed 
bonding loops is at the basis of both theoretical approaches. The W prediction for can 
thus be consistently used in connection with the FS approach [21] to predict the T and p 
dependence of the cluster size distributions. In the present case, the number of clusters 



4 



(per unit volume) containing / bi-functional particles and n three-funtional ones can be 
written [21] as 

Pnl = P3 —[P3Pbi^ - Pb)r[P2PbfWnl (3) 

PsPb 

_ (l + 3n-n)\ 
"^"^ ~ /!n!(n + 2)! 

where = 3N^/{2N2 + 3N^) and P2 = ^ — Ps are the probabilities that a randomly chosen 
site belongs to a three-functional or to a bi-functional particle, ph is given by Eq. [2| and 
Wni is a combinatorial contribution [21]. Distributions are normalized in such a way that 
Sjni+n>o(^ + '^)P"' = p2+p3- As showu in Fig. |3j-(top), on decreasing T, the pni distribution 
becomes wider and wider and develops a power-law tail with exponent -2.5, characteristic 
of loop-less percolation [22j. On further decreasing T, the distribution of finite size clusters 
progressively shrinks, since most of the particles attach themselves to the infinite cluster. 
Data show that Eq. |3| with no fitting parameters, predicts extremely well the numerical 
distributions at all state points, both above and below percolation. 

The three-functional particles act as branching points of the network formed by long 
chains of two-functional particles. Visualizing the structure of the system in term of chains of 
two-coordinated particles providing a link between the three-coordinates ones, it is possible 
to predict the number of finite size clusters composed of n three-functional units, irrespective 
of the number of bifunctional units. The system can thus be considered as a one-component 
fluid of three-functional particles forming clusters, in which the bonding distance between 
the three-functional particles is given by the length of the chains formed by the bi-functional 
units. Following again FS , it is possible to predict the ph value at which the systems develops 
a percolating structure: when Pb > pI = = 0.9256, an infinite cluster is present in the 
system. The percolation line is thus the locus of points in the phase diagram such that 
Pb{T,p) = pI, with pi,{T,p) given in Eq|2] Along the percolation locus, the product pA 
{pe^^^ at low T) is constant (from Eq|2|. 

It is also possible to predict the length distribution of chains of bi-functional particles. The 
number of chains (per unit volume) of length /, normalized in such a way that J^'i^i ^Pi — P2, 
is 

Pi = P2il - P2Pb)\p2Pby~^ (4) 

The pi distribution is thus always exponential. At low T, when ph ^ 1, the distribution 
becomes controlled only by the relative fraction of two-functional particles, providing a 



method for tuning the porosity of the fully connected percolating structure via the relative 
composition of the binary mixture. Indeed, when p;, — > 1 , the average distance / between 
branching points in the network becomes only a function of p2, i.e. I = 1/(1 — ^2)- Fig- ^ 
(bottom) shows that also the pi distribution is perfectly described by the combined W-FS 
theory. 

In the framework of FS approach it is also possible to evaluate the number density of 
clusters pc = J2ni Pm ^ function of pb, irrespectively of the cluster size. Below percolation, 
in the absence of bonding loops, the relation between p^ and pt, is linear, since each added 
bond decreases the number of clusters by one. Above percolation the relation crosses to 
a non-linear behavior, so that the number of clusters becomes one when pf, = 1. Within 
the FS theory, pc can be calculate for all pb values, assuming that finite clusters do not 
contain closed loops. As shown in Fig. |4| the simulation data conform perfectly to the 
theoretical expectation both below and above percolation. This suggests that, when the 
average functionality is small and the chains of bi-functional particles have a large persistence 
length, bonding loops in finite size clusters can be neglected. This agreement, which covers 
the entire range of pb values, implies that closed loops of bonds are statistically less favored 
than the corresponding open structure. The relative statistical weight results from the 
competition between the energy of forming the extra additional bond and the reduction of 
entropy associated to the closure of the loop. Hence, we interpret the absence of closed loops 
as resulting from the large configurational entropy of the long bi-functional chains. 

To further check the quality of the theory, we numerically evaluate the connectivity prop- 
erties of each studied state point, searching for the presence of clusters which are infinite 
under periodic boundary conditions. A state point is considered percolating when, account- 
ing for periodic boundary conditions, an infinite cluster is present in more than 50% config- 
urations. The resulting partitioning of the state points into percolating and non-percolating 
ones (see Fig. |5]) is fully consistent with the theoretical prediction of the percolation line 
Pb{T,p) = pI. The theory also predict a line of constant volume specific heat Cy maxima 
(also observed in chain polymerization [HI [23]), provided by the inflection point in the pb 
vs. T curves (Fig. [2]), which also agrees very well with the simulation results. The line of 
Cy extrema is also shown in Fig. |5j The presence of a maximum in Cy is a characteristic 
of bond-driven assembly and the locus of maxima in the T-p plane is one of the precursors 
of the self-assembly process for low functionality particles. 
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The W theory predicts a hquid-gas phase separation at small p for any non-vanishing 
amount of branching point [2H ES]- According to the theory, at low T and sufficiently small 
p the system phase separates into two phases of different density and connectivity. The 
theoretical spinodal curve, the line separating the stable (or metastable) state points from 
the unstable ones, is defined as the locus of points such that the volume V derivative of 
the pressure P vanishes, i.e. {dP/dV)T = 0. For the present model, it is located below the 
percolation line and the two lines merge asymptotically for T ^ and p ^ as shown in 
Fig. [5j The analysis of the numerical configurations for the two investigated state points 
which happened to be located inside the spinodal are indeed characterized by a bimodal 
distribution of the density fiuctuations and a very large value of the small-angle structure 
factor, indicating a phase-separated structure. This confirms that, in the absence of bonded 
rings, also the region of liquid-gas instability can be predicted using the W theory. The role 
of bond-rings on the thermodynamic of thermoreversible gels has been recently discussed in 
Ref. [ig. 

To summarize, this work provides the first fully solvable example of an equilibrium self- 
assembly process which goes beyond the equilibrium chain polymerization, complementing 
recent numerical studies and attempts at direct comparison between theory and simulation 
in related systems [151 ISSl |27]. The presence of a non negligible number of three- valent 
particles brings in a percolation phenomenon and a gas-liquid instability. These two loci 
are located in the region of T below the Cy maximum, suggesting a cascade of phenom- 
ena characterizing the self-assembly process in systems with low (M): onset of bonding 
(Cy max), percolation and, eventually at low p, phase separation. The resulting phase di- 
agram bears strong similarities with the ones discussed in previous mean field [28l [29] and 
numerical |15l |26l [30] studies of thermoreversible gelation. 

The possibility of theoretically describing the self-assembly process in the present model, 
designed to minimize the closed-bond loops effects, provides a benchmark for testing novel 
approaches and approximations accounting for more general cases. The quality of the pre- 
dictions makes it possible to exploit a fine tuning of the percolation properties. Using the 
recently synthesized colloidal particles [311 E2] or functionalizing colloids [331 Ell |35j it will 
be possible, by choosing the appropriate p and T, to generate materials with specific connec- 
tivity lengths, gels with desired porosity or build on the fractal nature of the self-assembled 
clusters. The material structure can be made permanent, if needed, by freezing the bonding 
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pattern with a fast cooling process. 

Finally, Pb appears to be the effective parameter controlling the structure of the system, 
since configurations at different p and T but equal pb are characterized by the same dis- 
tribution of cluster sizes. Configurations only differ in the relative distance between these 
clusters and, more interesting, on the lifetime of the resulting bonds, r ~ g^o/fesT^ 
thus in principle possible to decouple the effects related to the persistence of the clusters 
as units and the dynamics of clusters. The possibility of interpolating in a continuous and 
structural preserving way from chemical to physical gels may offer a unifying picture of these 
two apparently different arrest processes. 
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Figure Captions 



FIG. 1: Representation of the studied system for three different temperatures (in units of uo) at 
density pa^ = 0.04. Equally colored particles belong to the same cluster. The pictures show the 
structure of the system below percolation (left, ph < pf), where particles are aggregated in small 
finite clusters, at percolation (center, pf, Ri pf), where a spanning cluster first appears, and well 
above percolation (right, pb 3> p^), where most all of the particles belong to the infinite cluster. The 
mesh size of the fully connected system — provided by the average length of the chains connecting 
the three-functional branching points — can be finely controlled by the relative composition of the 
mixture. 

FIG. 2: Temperature (in units of uq) and density dependence of the bond probability pb- Points are 
simulation results based on Monte Carlo simulations for nine different densities and eight different 
temperatures. Lines are parameter-free predictions based on the Wertheim theory. Note that at 
low T the system reaches a fully bonded configuration. 

FIG. 3: Cluster size distribution for pa^ = 0.04 for some of the investigated T. Points are simulation 
data and lines are the corresponding theoretical curves, (top) The number of finite clusters (per 
unit volume) containing I bi-functional units plus n three-functional ones, (bottom) The number of 
chains of length I (per unit volume) composed of bi-bunctional particles. At low T, the distribution 
approaches its limiting (concentration controlled) value pi = P2P3P2~^, shown by the dashed line. 
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FIG. 4: Relation between the number of finite-size clusters Nc (irrespective of their size) and 
the bond probability. Symbols: simulation results for all studied densities. Solid line: Flory- 
Stockmayer predictions (both below and above percolation). Below percolation, in the absence 
of bonding loops, Nc is given by the difference between the total number of simulated particles 
N = 6000 and the number of bonds Nj,, since each added bond decreases Nc by one, i.e. A^^^ = 
N — Nh = N — ph^^^^^^^^ . This linear relation (dashed line) can be extended above percolation 
(dotted line), but never beyond the point where Nc < 1. The inset enlarges the region of large pb 
values, to provide evidence that the FS approach is valid over the entire pj, range. 



FIG. 5: Phase-diagram of the studied model. Lines are theoretical predictions: the solid (green) 
line is the spinodal curve obtained from the Wertheim equation of state, i.e. {dP/dV)T = 0; 
the dashed (magenta) line is the locus of points such that the specific heat is maximum, i.e. 
{dCv/dT)v = 0; the dotted (red) line is the percolation locus ph{p,T) = p'^ . Points are simula- 
tion results: the (black) triangles are the non-percolating equilibrium state points and the (red) 
squares are the equilibrium percolating state points. Gircles are (non-equilibrium) state points 
characterized by a bimodal distribution of the density fluctuations, indicating a phase-separated 
structure. 
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